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ABSTRACT 

We study the relation between nuclear massive black holes and their host spheroid gravitational 
potential. Using simple models, we analyze how gas is expected to be transported in the nuclear regions 
of galaxies. When we couple it with the expected gas lifetime given by the Kennicutt-Schmidt Law, this 
naturally leads to the 'Mbh - Mvidai' and 'Mbh - cr' relations. 
Subject headings: quasars: general - galaxies: formation - black hole physics 



L INTRODUCTION 

Most nearby massive spheroids (elliptical and spiral 
bulges) host nuclear massive black holes (MBH), whose 
mass scales as Mbh oc a'^~^ (the 'Mbh - cr' relation; Fer- 
rarese & Merritt 2000; Gebhardt et al. 2000). In recent 
years, several theories have been proposed to clarify the 
origin of this relation (Silk & Rees 1998; Blandford 1999; 
Fabian 1999; Burkert & Silk 2001; Zhao, Haehnelt & Rees 
2002; King 2003; Adams et al. 2003; Miralda-Escaude & 
KoUmeier 2005; Sazonov et al. 2005; Begelman & Nath 
2005). However, very little effort has been devoted to 
study how angular momentum is transfered in the cen- 
ters of galaxies. This is an unavoidable challenge for the 
understanding of the formation of MBHs and the origin of 
the 'Mbh - cr' relation. 

The study of MBH growth by gas accretion is usually 
focused on the study of accretion disks. However, these ac- 
cretion disks are Keplerian by nature and therefore have 
neglible masses compared to that if the MBH. They must 
be continously replenished, otherwise the mass of the MBH 
will not have a considerable growth. The key question in 
the growth of MBHs by accretion, is how to remove the 
large angular momentum of gas in a galaxy in order to 
funnel it into the accretion disk in the central sub-pc re- 
gion: the so called 'Fueling Problem'. Without a proper 
understanding of this 'Fueling Problem', specially in the 
last few 100 pc, it is impossible to understand the growth 
of MBHs and to elucidate why their masses correlate with 
the host spheroid properties. 

2. MASS TRANSPORT IN GALACTIC NUCLEUS 

Several mechanisms for fueling gas from R ^10 kpc 
down to 100 pc have been proposed in the literature, in- 
cluding galaxy interactions, stellar bars and gravitational 
resonances (see Shlosman, Begelman & Frank 1990 for a 
classic review, or Wada 2004 for a more recent one). How- 
ever, gravitational torques in galaxy mergers arise as the 
dominant process for fueling large amounts of gas down to 
the central few hundred parsecs and for triggering most of 
the MBH growth. In a merger, after a violently relaxed 
core is formed at the center, most of the gas will settle in 
a nuclear disk (with typically several hundred parsecs in 
diameter) that is rotationally supported against the over- 
all gravitational potential (Barnes 2002). The MBHs will 



migrate to the center, and merge in a timescale relatively 
short compared to the lifetime of the nuclear disk (Escala 
et al. 2004, 2005). 

2.1. Turbulent Disk 

Nuclear gas disks are characterized by being self- 
gravitating. It has been debated whether gravitational 
instability leads to disk fragmentation ('disk of clouds'; 
Shlosman et al. 1990) or to transient density fluctua- 
tions ('turbulent disk'; Wada 2004), depending on cool- 
ing/heating rates and thermal physics. I will mainly focus 
in the latter case because it is possible to advance further 
with analytical estimates. 

A turbulence-dominated disk can be modeled on a zero- 
order approximation as a steady a-disk. We choose to pa- 
rameterize the viscosity v in terms of the turbulent speed 
vturb instead of the sound speed, because the turbulence in 
a nuclear disk is expected to be highly supersonic (Downes 
& Solomon 1998) and therefore vturb becomes dynamically 
more relevant. In such a case, the disk height h is vturb/^ 
(being £7 the angular velocity) and the viscosity v, pa- 
rameterized by Shakura and Sunyaev as afl/i^, becomes 
V = av^^,^f^/il. Therefore in a turbulence-dominated a- 
disk with a flat rotation curve, the mass accretion rate 
can be written as 
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where Q = vturb^/TrGSg is the 'turbulent' Toomre's pa- 
rameter and Eg is the gas surface density. Eq. (1) repre- 
sents the mass accretion rate onto the inner accretion disk 
around the MBH (fueling). 

In nuclear disks the turbulence can be driven by two 
mechanisms: feedback from star formation and self- 
gravity. The thickness of these nuclear disks is indepen- 
dent of whether or not the disk hosts a starburst ( Jogee et 
al. 2005), suggesting that the turbulence is driven by self- 
gravity. We will now focus on the case where the turbu- 
lence is driven by gravitational instabilities and the energy 
is extracted from the disk rotation. Unfortunately, there is 
no consensus on the energy spectrum E(k) ex k"^ for com- 
pressible supersonic turbulence (Ballesteros-Paredes et. al 
2006). Boldyrev (2002) has presented a theoretical model 
suggesting that n « —1.74 and that appears to be con- 
firmed numerically ( Joung & Mac Low 2006 and references 
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therein). Boldyrev's theory suggests that in the inertial 
range, compressible supersonic turbulence obey the Kol- 
mogorov Law. In such a case the energy is transported 
in a cascade from the largest unstable scale characterized 
by the disk radius Rd and rotation velocity Vj-ot, to the 
turbulent eddies' scale characterized by the size Aturb and 
velocity vturb- This can be expressed as 
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We shall also consider that the gas is in a disk in rota- 
tional support against the overall gravitational potential. 
For simplicity we assume that v^ot and a are constant with 
radius, that is only strictly valid for an isothermal sphere. 
On the central kiloparsec of galaxies, stellar bulges are 
the dominant components of the gravitational potential 
even in the most massive gaseous nuclear disks known so 
far, where observations (Downes & Solomon 1998) found 
typically 5 times more mass in the bulge that in the gas 
(Mgas/Mstar ~ 1/5)- Evcn in such an extreme case, the 
rotation velocity is Vj-ot = (1-2 GMgtar/R)^^^ = 1.09 ct « a 

Finally, combining Vj-ot ~ cr with Eqs. (1) and (2) we 
find that the mass accretion rate onto the inner Keplcrian 
accretion disk is 
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Equation (3) shows that the mass transport in the nu- 
clear regions of galaxies is controlled by the spheroid veloc- 
ity dispersion a (the other factors are either dimensionlcss 
or physical constants). Therefore the rate at which the in- 
ner accretion disk around the MBH is replenished, strongly 
depends on the spheroid that hosts the MBH. 

2.2. Disk of Clouds 

A disk composed by randomly moving clouds is the sec- 
ond possible model for the nuclear disk. The accretion 
flow in this 'disk', can be expressed in the same form 
of Eq. 1 replacing vturb by the cloud velocity dispersion 
Av (Goldreich & Tremaine 1978). The formation of the 
clouds is driven by gravitational instabilities and their or- 
bital energy is extracted from the disk rotation, as in the 
gravity-driven turbulence case. Therefore we expect that 
Av scales with Vjot in a similar way as vturb does, leading 
to an expression for the mass accretion rate analogous to 
Eq. 3. However, the exact value of Av is highly uncertain 
and numerical simulations are required for an accurate es- 
timation. 

3. CLOSURE EQUATION: GAS LIFETIME 

In order to determine the total mass accreted, we need 
to know for how long gas will be available to be accreted 
by the MBH. 

Usually, AGN feedback is considered as the most likely 
physical process that determines the characteristic gas life- 
time, expelling most of the gas in the host galaxy when 
the Eddington Limit is exceeded. However, in an inho- 
mogeneous medium this limit can be greatly exceeded 
because photons tend to diffuse trough the low density 
medium. Several instabilities will produce inhomogeneities 



large enough to guarantee that the Eddington Limit is 
no longer valid, such as photon-bubble (Begelman 2001) 
or Kelvin-Hemboltz instabilities (Krumholz et. al 2005). 
Moreover on galactic scales, the ISM is endeed highly in- 
homogeneous and is hard to believe that most of the gas 
in a galaxy will be expelled due to AGN feedback. 

On the other hand, these massive nuclear disks are char- 
acterized by an ongoing starburst and therefore gas deple- 
tion due to star formation will eventually be relevant in 
the absence of other processes. Although star formation 
is still poorly understood, it can be parameterized by a 
well constrained empirical relation called the Kennicutt- 
Schmidt Law (Kennicutt 1998) which can be stated as 



Sgas = 0.017 



'^dyn 



(4) 



where Egas is the star formation rate per unit area, Egas is 
the gas surface density and Tdyn = Rd/vrot is the dynami- 
cal (i.e. orbital) timescale of the gas. These nuclear disks 
have typical radius of Rd ^ 0.1 Re, being Re the effective 
radius of the host spheroid. As shown in §2, the bulge 
dominates the overall gravitational potential (v^ot ~ o") 
and therefore the dynamical timescale of the gas can be 
expressed as Tdyn = 0.1 Ro/c. Thus Eq. 4 implies that the 

gas lifetime tgas, defined as Sgas/Egas, is tgas = ^Re/c^ 
where ^ « 5.9. 

The final mass of the MBH will be determined by 
Mbh — rjM tgas , where rj is the fraction of gas in the inner 
accretion disk that finally ends inside the event horizon of 
the MBH. Using Eq. (3) the final MBH mass becomes 
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which reproduces the correlation between Mbh and the 
bulge mass Mviriai = kRetr^/G (Marconi & Hunt 2003; 
Haring & Rix 2004), being k=8/3 for an isothermal sphere. 

Using the Faber- Jackson relation Re oc a'^ in Eq. 5, we 
obtain the standard 'Mbh - f' form (Mbh oc cr^). Since 
the 'Mbh - cr' and 'Mbh - bulge virial mass' correlations 
have similar levels of scatter, the Faber-Jackson relation 
should not have larger scatter than either one in the galax- 
ies where the correlations are measured. 

Finally, we need to evaluate the zero point in Eq. 5. For 
gravity-driven turbulence a may be as large as 0.3 (Gam- 
mie 2001). In nuclear disks vturb/vrot ^ 1/5 (Downes 
& Solomon 1998), implying that Aturb/Rd ^ 0.008 (Eq. 
2). Using a typical value of Q = 1.5 (Martin & Kenni- 
cutt 2001), we find Mbh ~ 0.009 ??Mviriai- To reproduce 
the 'Mbh - bulge virial mass' relation, we need to have 
7] w 1/6. Thus, we need that less than 20% of the mass 
fed onto the inner accretion disk ends up inside the event 
horizon of the MBH. 

4. NUMERICAL STUDY 

The conclusions found in §2 and §3 are based in a steady- 
state model (a-disk). However, previous numerical stud- 
ies on the mass transport in the inner 100 pc of galaxies 
(Wada & Norman 2002; Wada 2004) found that the mass 
accretion rate is highly time dependent with fluctuations 
of several order of magnitudes. For that reason, we have 
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designed a simple numerical experiment to test the validity 
of the a-disk model as a zero-order approximation. 

In massive nuclear disks, self-gravity, cooling and feed- 
back effects from star formation play a major role in de- 
termining the physics of the ISM and therefore must be 
included to accurately model the angular momentum and 
mass transport in the disk. However, a fully realistic model 
has multiple relevant parameters to be studied (such as the 
ones in the star formation prescription, heating/cooling 
processes, etc) and requires a long exploratory numerical 
study. This is beyond the scope of this Letter and will be 
studied in subsequent papers. For that reason, we decided 
instead to adopt an idealize model where star formation 
is not included and the thermodynamics is simplified by 
using an adiabatic Equation of State (E.O.S.), but where 
the velocity dispersion of the bulge a is the only varying 
parameter. With this E.O.S. we will not be able to com- 
pute accurate mass accretion rates, but we will be able 
to compare the changes in the accretion rates for differ- 
ent values of a and then reveal the role of a in the mass 
transport due to gravity-driven turbulence. 

The model consist in a Mg = 5 x 1O^M0 gaseous disk 
around a MBH (Mbh = 10%©) that is embedded in a 
stellar bulge. The MBH and the bulge are modeled as 
time-independent external potentials. For the bulge we 
use the potential of an isothermal sphere for r > 12pc, 
and the potential of an homogenous sphere for r < 12pc. 
We modeled the MBH by a Plummer potential with core 
radius of 1 pc. We solved the hydrodynamic and Poisson 
equations using the FLASH adaptative mesh refinement 
hydrodynamics code (Fryxell et. al. 2000). The carte- 
sian grid covers a 256^ x 32 pc'^ region around the galactic 
center, with a spatial resolution of 2pc. We use outflow 
boundary conditions. 

The initial condition is a rotationally supported disk 
with uniform density profile, radius of 125 pc and thick- 
ness of 4 pc. The initial temperature is set to lO'^K over 
the whole region. Random density and temperature fluc- 
tuations of less than 1% are added to the initial uniform 
disk. We performed four different runs where we vary the 
mass of the bulge, in such a way that the velocity disper- 
sion of the bulge a is: 100, 200, 300, and 400 kms"!. For 
a bulge modeled by an isothermal sphere, our rotationally 
supported disks have flat rotation curves. For example the 
run with the highest gas fraction Mgas/Mgtar has a rota- 
tion velocity Vj-ot of 109 km s~^. This corresponds to the 
a = 100km s~^ run and the Mgas/Mgtar ratio is similar 
to those found in the most massive gaseous nuclear disks 
(Downes & Solomon 1998). 

4.1. Results 

Figure 1 illustrates a representative stage in the evolu- 
tion of the system, showing the face-on density distribu- 
tion at the plane of the disk for run 100, at the time t = 
6 Myr. The flgure shows a complicated multiphase and 
highly turbulent structure, in qualitative agreement with 
the findings of Wada & Norman (2002). 

As mentioned before, we start our simulation with an 
uniform disk with random fluctuations of less than 1% 
which are non-linearly amplified by self-gravity. The sys- 
tem is then let to evolve approximately for an orbital time 
torb = 2Rd/vi.ot7 in order to fully develop the turbulence 
in the disk. Because torb varies among the different runs. 



we choose the longer torb as starting point for our com- 
parison. This corresponds to the run with a = 100km s^^ 
which has an orbital time torb = 2.44Myr. 

Figure 2a shows the time evolution of the gas mass inside 
a radius of 2 pc around the MBH (our resolution limit), 
for the four different values of a: 100 (black), 200 (red), 
300 (blue) and 400 kms^^ (green). This figure shows a 
strong dependence on the bulge velocity dispersion, and 
the final gas mass (r < 2pc) varies by almost a factor 70 
among the runs. In all the runs, the gas mass (r < 2pc) 
shows an approximately exponential growth. 

In Fig. 2b, we plot the time evolution of the mass ac- 
cretion rate inside 2 pc, M(r<2pc), for the different runs 
and using the same color code as in Fig. 2a. The accre- 
tion rates are highly variable on timescales of lO^yr and 
with fluctuations of a couple of orders of magnitude in M, 
as found by Wada (2004). However, we found accretion 
rates considerably lower than those of Wada (2004) . The 
origin of this discrepancy rely in the adiabatic E.O.S. This 
assumption produces unrealistically high Toomre's Q pa- 
rameters compared to the observed ones (Q ~ 1.5; Martin 
& Kennicutt 2001), decreasing the mass accretion rates on 
Eq. 3. 

Finally, in Fig. 3 we plot the average mass accretion rate 
< M > as a function of a. We compute < M > simply 
as [M(tf) — M(ti)]/[tf — ti], where ti and tf are the initial 
and final times. We choose ti = 2.4 Myr to ensure that 
all the runs have evolved for at least an orbital time torb 
and at that time all the runs have almost the same gas 
mass inside 2 pc. The final time of the runs is tf = 7 Myr. 
It is found that the average mass accretion rate < M > 
is strongly dependent on a, as seen in Fig. 3. The black 
line is a least squares fit to the points in Fig. 3 and cor- 
responds to < M >cx (T^-^. Therefore < M > scales with 
a as expected from Eq. 3, despite the complex behavior 
of the gas (Fig. 1) and that the mass accretion rates are 
highly time-dependent (Fig. 2). 

5. DISCUSSION 

We perform a simple numerical experiment that illus- 
trates the role of the spheroid gravitational potential on 
the mass transport in a gravity-driven turbulent disk. We 
find that -on average- the mass accretion rate behaves as 
predicted by the zero-order approximation (M ex ct^), how- 
ever our adiabatic simulations underestimate the mass ac- 
cretion rates. 

Using more realistic models by Wada & Norman 
(2002) that found M - 1.4M0yr-i for an spheroid 
with (T — lOOkms"^, we estimate that M « 
11.2(cr/2OOkms"i)3M0yr~i (Eq. 3). If we assume that 
10% of the mass fueled into the inner accretion disk is 
accreted by the MBH (77 = 0.1) and use the 'Mbh - 
cr' relation, we find that MBHs could gained a mass 
Mbh in a rapid growth phase of tgrowth = Mbh/'?M « 
10^ (Mbh/10* Mq)^/^ yr. This is in agreement with recent 
observations (Borys et al. 2005), that are consistent with 
a model where MBHs in submillimeter galaxies undergo a 
rapid growth to reach the local 'Mbh - cr' relation. 

We find that the 'Mbh - cr' relation arises naturally from 
the fueling of either a turbulent disk or a 'disk of clouds' 
(coupled with the Kennicutt-Schmidt Law), and therefore 
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is hard to distinguish between the two models. However, 
the zero point could be different allowing to discriminate 
between both disks. Also, the zero point allow us to ad- 
dress the possible role of star formation driven turbulence. 
Detailed numerical simulations of these various disks are 
necessary to determine a realistic zero point that could be 
compared with the one measured in the 'Mbh - f' relation. 
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Fig. 1. — Density distribution in the plane of the gas disk, coded on a logarithmic scale, at time t = 6Myr for the run 100. The black dot 
indicates the position of the MBH. The medium is characterized by high density clumps and filaments, embedded in a less dense medium. 
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Fig. 2. — (a) Time evolution of the gas mass inside 2 pc around the MBH (our resolution limit), for the four different values of cr: 100 
(black), 200 (red), 300 (blue) and 400 kms~^ (green), (b) Time evolution of the mass accretion rate inside 2 pc, M(r<2pc), using the same 
color code as in (a). 
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Fig. 3. — Mass accretion rate averaged over the time interval [2.4Myr,7Myr], for the four different values of cr: 100, 200, 300 and 400 
kms^^ (in the same color code as in Fig. 2). The black line is the best-fitting of the points using minimum least squares, and corresponds 
to < M XX cr^-i. 



